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A MASSIVE PROTOSTAR FORMING BY ORDERED COLLAPSE OF A DENSE, MASSIVE CORE 
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ABSTRACT 

We present 30 and 40 imaging of the massive protostar G35.20-0.74 with SOFIA-FORCAST. 
The high surface density of the natal core around the protostar leads to high extinction, even at 
these relatively long wavelengths, causing the observed flux to be dominated by that emerging from 
the near-facing outflow cavity. However, emission from the far-facing cavity is still clearly detected. 
We combine these results with fluxes from the near-infrared to mm to construct a spectral energy 
distribution (SED). For isotropic emission the bolometric luminosity would be 3.3 x 10 4 L . We 
perform radiative transfer modeling of a protostar forming by ordered, symmetric collapse from a 
massive core bounded by a clump with high mass surface density, E c i. To fit the SED requires 
protostellar masses ~ 20 — 34 M depending on the outflow cavity opening angle (35° — 50°), and 
E c i ~ 0.4 — 1 g cm -2 . After accounting for the foreground extinction and the flashlight effect, the true 
bolometric luminosity is ~ (0.7 — 2.2) x 10 5 L©. One of these models also has excellent agreement 
with the observed intensity profiles along the outflow axis at 10, 18, 31 and 37 /mi. Overall our results 
support a model of massive star formation involving the relatively ordered, symmetric collapse of a 
massive, dense core and the launching bipolar outflows that clear low density cavities. Thus a unified 
model may apply for the formation of both low and high mass stars. 
Subject headings: stars: formation 



1. INTRODUCTION 

Massive stars impact many areas of astrophysics. In 
most galactic environments they dominate the radiative, 
mechanical and chemical feedback on the interstellar 
medium, thus regulating the evolution of galaxies. Many 
low-mass stars form in clusters along with massive stars, 
and their protoplanetary disks can be affected by this 
feedback also. There is some evidenc e that our ow n so- 
lar system was influenced in this way (Adams 2010). De- 
spite this importance, there is no consensus on the basic 
formation mechanism of massive stars. Theories range 
from Core Accretion, i.e., a scaled- up version of low- 
mass star formation (e.g. the Turbulent Core Model of 
iMcKee fc Tanll2003l . MT03), to Co mpetitive Accretion at 
the centers of form ing star clusters ([Bonnell et al.|[200Tl 
Wan g et al.l 12010). to Stellar Collisions (jBonnell et al.l 
1998). 

The Core Accretion theory predicts the existence of 
an envelope-fed, central accretion disk and relatively or- 
dered and collimated bipolar outflows powered by accre- 
tion around a massive protostar. The other formation 
mechanisms predict a much less ordered geometry of gas 
and dust surrounding the protostar. 

Collimated bipolar ou tflows have been obs erved from 
massive protostars (e.g. iBeuther et al.l [20021 ). Outflows 
may limit the formation efficiency (mass ratio of final 
star to i nitial gravitationally-bou nded protostellar core) 
to ~0.5 (Ma tzner fc McKea [T999), since they expel ma- 



terial along polar directions. The resulting low-density 
cavities have been proposed to significa ntly affect the ap- 
pearance of the proto star in the MIR (|De Buizerl 2006, 
iDe Buizer et al.l 12012), and this is seen in the ra diative 
tran sfer (RT) cal c ulatio ns of IZhang fc fanl ([201li ZT11) 
and IZhang et all (12011 ZT M13) using the RT code de- 
veloped bv lWhitnev et all (|2003f ). 

The Turbulent Core Model also relates core properties 
to the surrounding self-gravitating star-cluster-forming 
clump. For a marginally unstable pre-stellar core of a 
given mass, its size, density and subsequent accretion 
rate all depend on the pressure set by the mean mass 
surface density of the clump , E c i. Observed valu es of E c i 
are-O.l-lgcm- 2 fMT03: lButler fc Tan|[2012h . leading 
to high predicted accretion rates, ~ 10 -4 — lO _3 M yr _1 , 
potentially important for overcoming the high radiation 
pressure from the massive protostar. 

G35.20-0.74 (hereafter G35.2) is a massive protostar 
associated with a well-defined outflow. At a distance 
of 2.2 kpc ([Zhang et al.l l2009h , radio continuum emis- 
sion shows that it contains an ultracompact HII region 
at the center in the f orm of a collimated jet along the 
north -south direction (He aton fc Little 1988, G ibb et all 
2003), which may be an examp le of an Outflow- Confined 
HII Region ([Tan fc Mc Kee 20 0l). Th e outflow is also ob- 
served in the NIR ([Fuller et al.ll200ll ) with the northern 
cavity brighter than the southern one. Ground-based 10 
and 18 jam observations only reveal the northern cavity 
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(|De Buizerl [2006), suggesti ng it is inclined towards us. 
Outflow is also seen in CO (jGibb et al.ll2003llBirks eTall 
2006), including a wider angle flow in the NE-SW direc- 
tion. It is unclear if this results from a separate driving 
source or is related to the N-S jet, but misaligned by in- 
teraction with the surrounding core/clump. In previous 
studies, an early B star wa s thought to be present at 
the center (|Dent et al.lll985l). surrounded by an envelope 
with mass ~ 500 ((Little et all 1 1991 to - 3800 M© 
(jParon fc Weidmannll2010D . 

Here we present the first SOFIA-FORCAST images of 
G35.2. Combined with other multiwavelength data, we 
construct the SED. We present a simple exploration of 
parameter space of RT models of massive protostars and 
use the observed SED and intensity profiles along the 
outflow axis to constrain the models. 

2. OBSERVATIONS & DATA REDUCTION 

The SOFIA-FORCAST (jHerter et al.l f2012h observa- 
tions of G35.2 were performed on May 24 and 26 (UT) 
2011 in the 11.3, 19.7, 31.5, and 37.1 /xm filters at an 
altitude of 43000 ft. The chopping secondary of SOFIA 
was driven at 4 Hz, with a matched chop and nod throw 
of 60", and with a nod performed every 30s. The final ef- 
fective on-source exposure times were 909 s at 11.3 /im, 
959 s at 19.7 /im, 4068 s at 31.5 /im, and 4801 s at 
37.1 /im. The fluxes were calibrated by the SOFIA data 
reduction pipeline, which during the Basic Science pe- 
riod h as an estimated 3a < 20% calibration error in all 
filters (jHerter et al.ll2012l ). 

The 11 and 31 /im images were observed simultane- 
ously with the dichroic mode of FORCAST; and simi- 
larly for the 19 and 37 /im images. The relative offsets 
between the short and long wavelength arrays are known 
to better than la = 0.15". The best registration of the 
19 /im im a ge wi th respect to the 18.3 /im image from 
iDe Buize r (2006) was obtained using a x 2 minimization 
algorithm. The latter was convolved with a Gaussian 
profile to reach the same resolution of the FORCAST 
image. Several features in the two images match up well 
after the registration increasing confidence in the tech- 
nique. We adopt the pixel size (0.768") as our conserva- 
tive la uncertainty in the relative offset between the 18.3 
and 19.7 jam images. Considering the absolute astromet- 
ric uncertainties of ~ 0.2" from the 18.3 /im image, we 
get an absolute astrometric uncertainty of 0.79" from the 
19.7 /im image. The absolute astrometric uncertainty of 
the 37.1 /im image is then estimated to be better than 
0.8". However, the 11.3 /im FORCAST image suffered 
from such low S/N that it could not be used to per- 
form similar registration for the 31.4 jam image. We had 
to use the x 2 minimization algorithm to find the best 
registration of the 31.4 /im image and the astrometri- 
cally calibrated 37.1 /im image, noticing that they look 
remarkably similar morphologically. The absolute astro- 
metric uncertainty in the 31.4 /im image is estimated to 
be 0.89". 

Using standard stars observed throughout several Ba- 
sic Science flights, an average FWHM for each wave- 
length was determined: 3.9" at 19.7 /im, 4.1" at 31.5 /im, 
and 4.5" at 37.1 jam. The images were then deconvolved 
using the m aximum likelihood method (jRichar dsonll 19721 
lLucvl[l974} ) , which yielded an approximately factor of two 
better resolution than the pre-deconvolved images. The 



deconvolved images compare favorably to simple unsharp 
masking of the original images, and in the case of the 
19.7 /am image, it compares well to the Gemini 18.3 jam 
image, and thus the substructures revealed in the decon- 
volved images are believed to be real. 

We also retrieved calibrated level 2 data (AOR 
1342241157, HOBYS GT program) from Herschel- 
SPIRE and a single long integration SCUBA 850/450 
jam jiggle- map form the JCMT archive at CADCQ. The 
SCUBA data wer e reduced and calibr ated in a standard 
way using SURF (jSandell et al.ll2001[ ). 

3. RESULTS 

Fig. (pQ) shows the FORCAST images at 31.5 and 
37.1 /am a long with the 11 and 18 /am images by Gemini- 
T-ReCS (De Buizer 2006) and 3.6 a nd 8 /am images 
by Sp itzer-IRAC (GLIMPSE Survey, IChurchwell etaD 
2009). The airborne FORCAST images have a dy- 
namic range of ~ 4 orders of magnitude, similar to the 
space-based IRAC observations and much higher than 
the ground-based T-ReCS observations. The origin of 
the images is chosen to be the radio source G35.2N 
(R.A.=18h58ml3.033s, DEC.=+ld40m36.14s). While 
the 11 and 18 jam images only show the northern side 
of the outflow, the FORCAST 31 and 37 jam images also 
reveal the fainter southern outflow cavity. This cavity is 
also seen in the IRAC images, again being fainter than 
the northern side. 

This trend agrees with the results of the RT sim- 
ulations by ZT11 and ZTM13, which showed that at 
~40 /im, the near-facing outflow cavity is the dominant 
feature but, because of lower extinction, the structure 
along the outflow axis is more symmetric than at shorter 
wavelengths. In these models, the far-facing outflow ap- 
pears in the NIR mainly due to the scattering by dust in 
and around the outflow cavity. 

With their high dynamic range, the FORCAST im- 
ages also show extended emission on larger scales, which 
is missing in the ground-base MIR images. The morphol- 
ogy of this extended emission has a resemblance to that 
seen with IRAC. This extended emission seems to be re- 
lated to the cavity wall and/or a wider- angle part of the 
outflow. It is also possible that some, probably small, 
fraction of this emission is produced by other luminosity 
sources. 

3.1. SED and Protostellar Model Global Parameters 

We integrated over a 40" x 50" region to find total 
31 and 37 /am fluxes of 527 and 922 Jy, respectively. No 
significant emission was seen outside of this area in ei- 
ther ban d. We constructed th e SED using data from 
2MASS (Skr utskie et~aTl l2006h. G LIMPSE. Gemini-T- 
ReCS, IRAS-IGA (jCao et al. 1997), Herschel- SPIRE and 
JCMT-SCUBA (this paper), and IRAM 30m Telescope 
(jMoonev et al.lfl995h . Total 2MASS and IRAC fluxes 
were integrated from background-subtracted images in 
the same region as the FORCAS T bands. The 1 and 
18 jam fluxes (3.06 and 46.87 Jy, IDe Buizerl l2QQ6h were 
integrated with a 30" aperture (again, there are no sig- 
nificant sources in the immediate vicinity). For IRAS 
bands, the beams are so large that the 40" x 50" re- 
gion does not contain all the flux from this source, so 

1 Guest user, Canadian Astronomy Data Center 
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Fig . 1. — Multiwavelength images of G35.2. From left to right, IRAC (GLIMPSE) at 3.6 and 8 /mi; T-Recs at 11 and 18 /mi (|De Buizerl 
2006 ); (deconvolved) F ORCAST at 31 and 37 /im. The white contours in the 11 /im image show 15 GHz radio continuum emission 
(Heaton Sz Little 1988). Total fluxes are labeled in the upper-right corners and maximum surface brightnesses in the lower-left. The dashed 
line indicates the axis of the outflow deduced from the MIR and radio emission. Circles indicate image resolutions. 



TABLE 1 

Parameters of the four fitted models. 



Model 


Core mass 

M c (M ) 


Mean surface density of 
the clump S c i (gem -2 ) 


Protostellar 

mass, ra* (Mq) 


Outflow opening 
angle 9 w ,esc 


Inclination 

■ a 

^view 


Foreground 
extinction Ay 


Luminosity 

^bol (£©)! 


Model 1 


240 


1 


34 


51° 


58° 


8 


2.2 x 10 5 


Model 2 


240 


0.7 


26 


45° 


51° 


8 


1.2 x 10 5 


Model 3 


240 


0.7 


22 


40° 


44° 


15 


9.0 x 10 4 


Model 4 


240 


0.4 


20 


35° 


43° 





6.6 x 10 4 



inclination between the line of sight and the outflow axis. 

b True bolometric luminosity after allowing for foreground extinction and flashlight effect. 



1Q -6| , , , , , 

- Model 1 

_ Model 2 




A( / u,m) 

Fig. 2. — SED of G35.2 (solid squares) with observed fluxes from 
2MASS H & K bands, IRAC from 3.5 to 8 /mi, T-ReCs 11, 18 /im, 
FORCAST 31, 37 /mi, IRAS 60, 100 /mi, SPIRE 250, 350 /mi, 
SCUBA 450, 850 /zm, and IRAM 1.3 mm. The lines show four 
model SEDs (see text). 



we use larger apertures (6.5' x 3') that encircle all the 
emission from the source while avoiding fluxes from sig- 
nificant nearby sources. We obtained fluxes at 60 and 
100 /im of 2.37 x 10 3 and 4.35 x 10 3 Jy. They may be 



considered as upper limits. For SPIRE and SCUBA im- 
ages, we have carefully subtracted the backgrounds with 
two or more component Gaussian fitting. The fluxes are 
546, 420, 123, and 16 Jy at 250, 350, 450, and 850 /im 
respectively with source sizes of 27" x 18" , 25" x 3.8" 
at 250 and 350 /im, and 22" x 12" in SCUBA bands. 
The SED is shown in Fig. ((2]). 

We use the SED to constrain the properties of the pro- 
tost ar, its core and the surrounding clump by fitting with 
the RT models of ZT11 and ZTM13. The models are de- 
veloped self-consistently for a core embedded in a high 
pressure clump forming a massive star via core accretion 
(MT03), including a treatment of rotating infall, an ac- 
tive accretion disk and a disk wind filling the outflow 
cavity. We have further extended the ZTM13 models to 
also include emission from the surrounding clump (see 
Fig. [3]). Since we have not yet constructed a full model 
grid, we limit ourselves to achieving only a moderately 
good fit by exploring the parameters manually and with a 
quite coarse sampling of parameter space. We show four 
models that have relatively good matches to the SED in 
Fig. (|2|) with their parameters listed in Tab. (pQ). 

Model 1 (Fig. [3]) contains a core of initial mass M c = 
240 Mq , and a protostar of mass m* = 34 M forming 
at its center. The core is embedded in a clump with 
E c i = 1 gem -2 , thus setting its radius R c = 0.11 pc 
(MT03). We distribute the clump material between 1 
and 10 R c with a power-law density profile (with index of 
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Fig. 3. — Cartoon showing geometry of Model 1. Density and temperature profiles are on the upper-left and lower-right sides of the 
outflow axis, respectively. White dotted lines show disk wind streamlines. The low temperature region close to the outflow axis is dust-free 
wind launched from the inner hot dust-free disk. The clump and the outflow extend to 10 core radii. The simulated images (60 /; x 60") of 
Model 1 are also shown. The scalebars are in \og(S/S avg ), where S avg is the average intensities of the near-facing outflow in each band (see 
Section [32]). For the T-ReCs and FORCAST bands, we apply same dynamic ranges as in Fig. {T}, while for the IRAC bands, we adopt 
larger dynamic ranges to show the far-facing outflow. 



-1.75, e.g. jMueller et al. 2002) that joins the core bound- 
ary smoothly and reaches E c i 1 g cm -2 by 10 R c 
(counting two sides along a central sight line). Then 
the original clump mass (assuming spherical symmetry, 
i.e. before some material is swept-up by the outflow) 
is 4.8 x 10 3 M Q . The size of the disk is set by as- 
suming a ratio of rotational to gravitational energy of 
2%, yielding an outer radius of 940 AU. Its mass is as- 
sumed to be 1/3 of the stellar mass. The accretion rate is 
5.1 x lO^Moyr -1 . The total luminosity is 2.2 x 10 5 L© 
including protostar, accretion boundary layer and disk. 
An outflow cavity, filled with a disk wind, including dust 
along streamlines launched from the cooler disk outside 
of 7 AU, is included. Following ZT11 and ZTM13, the 



cavity opening angle, 6 U 



51°, to yield a final star 



formation efficiency from the core of 50%. The view- 
ing angle between the sight line and the outflow axis is 
#view = 58°. At this viewing angle, if the received flux 
was assumed to be isotropic then a bolometric luminos- 
ity of 3.6 x 10 4 Lq would be implied. After considering 
the foreground extinction as a free parameter to better fit 
the NIR data (Ay = 8 mag for Model 1) and applying an 
appropriate aperture size (~ 50"), the inferred isotropic 
bolometric luminosity based on the received flux becomes 
3.3 x 10 4 Lq. Three other models with smaller values 
of E c i, ra* and # v iew are also shown in Fig. (|2j) (see also 
Table [Q). 

Overall, good agreement is achieved between these 
model SEDs and the observations. Small discrepancies 
in the NIR and IRAC bands may be due to patchy ex- 



tinction and/or poorly- modeled polycyclic aromatic hy- 
drocarbon (PAH) features and/or thermal emission from 
transiently heated small dust grains. At 60 fim the ob- 
served IRAS flux is lower than the model (after con- 
volving with the relatively broad filter response of IRAS 
60 jam band) by ~25%. These results show the SED fit- 
ting has a certain degeneracy in model parameter space, 
especially when the opening angle is not well constrained 
(ZTM13 found the outflow, especially the bright part, at 
10 to 20 fim appears narrower than the actual cavity, 
so the MIR morphology may not be a reliable tracer of 
the total outflow opening angle). Models with narrower 
outflow cavities reprocess larger fractions of their lumi- 
nosities and thus have smaller bolometric flashlight effect 
corrections for the considered viewing angles. Thus lower 
luminosity protostars are able to reproduce the observed 
FIR peak of the SED. Then, lower values of E c i and # v i ew 
are needed to match the MIR fluxes. Observations at 
10 /im with higher spectral resolution will provide strong 
diagnostics for testing these models. However, even with 
a large uncertainty in the outflow opening angle, the pro- 
tostellar mass is still estimated to be ~ 20 — 34 M and 
E c i ~ 0.4 — 1 gem -2 , i.e., G35.2 is a massive protostar 
forming from a massive core embedded in a high surface 
density clump. 

3.2. Resolved Intensity Profiles Along the Outflow Axis 

The images of Model 1 are shown in Fig. (|3j) for refer- 
ence. The other three models selected from SED fitting 
show slightly narrower outflow cavities. We compare the 
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Fig. 4. — Intensity profiles along the outflow axis. The squares 
are observational data sampled at intervals of the resolutions of 
the instruments (intervals of 2 Sz 3 x resolution are used for 18 
& 10 fim) with errors composed of systematic flux uncertainties 
(assumed to be 20%) and estimated background noise. The lines 
are model profiles. 



intensity distribution along the outflow axis predicted by 
these models with observations in Fig. (gj). The axis di- 
rection of G35.2 is chosen via the radio continuum and 
MIR morphology (dashed line in Fig. (pQ) with P.A. of 
6°). The model profiles are all convolved with the cor- 
responding instrument beams. At each offset from the 
center, we average over a perpendicular width of 2" to 
estimate profiles for both model and data. All model and 
observed profiles are normalized to the average values 
over the near-facing regions that have significant emis- 
sion. We also add constant background ambient intensi- 
ties, which may be either due to instrumental noise (i.e. 
in the T-ReCS data) or from additional ambient inter- 
stellar material. Note we have not attempted to use the 
lateral width of the cavities to constrain the models 

We emphasize that a detailed model search to fit these 
profiles separately from the SED has not been performed, 
but still Model 1 agrees very well with the observations 
at 37, 31, 18 and 10 /im, producing the right peak po- 
sitions and asymmetries of the two sides of the out- 
flow. At 10 /im, Model 1 predicts higher extinction to- 
wards the center. A possible explanation could be the 
strength (depth and width) of the silicate absorption fea- 
ture adopted in the dust model used in the RT calcula- 
tions (ZTM13). In the IRAC bands, while a reasonable 
fit is achieved for the near-facing outflow cavity, the in- 



tensities of the far-facing side predicted by the model 
are too low compared to those observed. Possible rea- 
sons for this include: (1) non-uniform, possibly patchy, 
extinction from the clump or the adopted uniform "fore- 
ground" extinction, that may actually be relatively local 
to the source (e.g. a foreground extinction of Ay = 15 
lowers 3.6 /im flux by a factor of two); (2) enhanced emis- 
sion from PAHs and transiently-heated small dust grains 
in and around the outflow cavity. Models 2, 3 and 4 also 
have good matches with the data in the near-facing and 
far-facing wings at 37, 31, 18 and 10 /im, but their pro- 
files become too peaked towards the center due to the 
lower extinction of their lower E cores and clumps. 

These results, especially the intensity profiles for Model 
1 from ~ 10 — 40 /im (for which the modeling has po- 
tentially fewer problems associated with small grain and 
PAH emission), again support the paradigm that G35.2 
is a massive protostar forming from a high E core and 
clump. The match of the SED and intensity profiles with 
observations becomes less good when E c i becomes lower. 
We estimate that E c i should be no lower than 0.4g cm -2 . 
The profiles also suggest that the accretion flow and out- 
flows are relatively well-ordered and symmetric. 

4. DISCUSSIONS & CONCLUSIONS 

SOFIA-FORCAST provides imaging with high dy- 
namic range similar to space-based instruments at a 
unique wavelength region ~30 to 40 /im, where lower 
extinction allows us to search for the predicted (ZT11, 
ZTM13) far-facing outflow cavity from a massive proto- 
star forming from a high surface density core. At longer 
wavelengths the emission is predicted to become even 
more symmetric, being dominated by cooler dust in the 
core/clump. Our SOFIA-FORCAST observations at 31 
and 37 /im did reveal emission from the far-facing out- 
flow cavity of G35.2, which was too faint to detect by 
ground-based T-ReCs 11 and 18 jam observations. 

We compiled the NIR to mm SED of G35.2. RT mod- 
eling of a massive protostar forming from a massive core 
bounded by a high E clump gave good agreement with 
this SED for four models. Depending on the outflow cav- 
ity opening angle (35° to 50°), we found m* ~ 22— 34M , 
L ho \ ~ (0.7 - 2.2) x 10 5 L© and E cl - 0.4 - 1 g cm" 2 . 
Model 1 also produced intensity profiles along the out- 
flow axis that fit the observations well at 10, 18, 31 and 
37 /im, without need for extensive fine tuning. These re- 
sults indicate G35.2 is a massive protostar, forming from 
high surface density core and clump, via relatively or- 
dered, symmetric collapse and accretion. Powerful bipo- 
lar outflows are being launched and have cleared wide- 
angle cavities, that are also relatively symmetric. 

A protostar with the luminosity estimated here (~ 
1 x 10 5 Lq) is expected to drive a CO outflow with 
mom entum flux of ~ 0.1 M0yr -1 kms -1 ([Richer et al.l 
2000), which is much larger than the observed value: 
even assuming that the larger NE-SW CO outflow is 
also driven by this source, the total momentum flux 
is still only ~ 003 M^yr^kms" 1 (jGibb et al.l 120031 . 
Bir ks et al.ll2006h . This may indicate that the outflow 
from G35.2N, being partly ionized as indicated by the 
observed radio continuum jet, may be relatively deficient 
in CO emission. Note that the modest misalignment of 
the larger CO outflow with the radio and MIR jet may 
result from the interaction of the wider angle part of the 
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outflow with a core/clump that does have some moderate 
lateral asymmetries. 

We conclude that, at least in the case of G35.2, mas- 
sive star formation appears to be well described by pre- 
dictions of the Turbulent Core Model (MT03), even with 
a highly idealized realization of the model utilizing per- 
fect axissymmetry and s mooth density structures (see 
iCunningham et all I2Q11I for a numerical simulation of 
this model including outflows). Such ordered accretion 
and outflow is not expected in models involving mas- 
sive star formation via competitive accretion or stellar 
collisions, nor indeed in cores that are highly turbulent. 
This may indicate that magnetic fields are dynamically 
important in setting the structure of the core. 
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